home *** CD-ROM | disk | FTP | other *** search
-
-
-
- SSSSGGGGEEEELLLLSSSSDDDD((((3333SSSS)))) SSSSGGGGEEEELLLLSSSSDDDD((((3333SSSS))))
-
-
-
- NNNNAAAAMMMMEEEE
- SGELSD - compute the minimum-norm solution to a real linear least squares
- problem
-
- SSSSYYYYNNNNOOOOPPPPSSSSIIIISSSS
- SUBROUTINE SGELSD( M, N, NRHS, A, LDA, B, LDB, S, RCOND, RANK, WORK,
- LWORK, IWORK, INFO )
-
- INTEGER INFO, LDA, LDB, LWORK, M, N, NRHS, RANK
-
- REAL RCOND
-
- INTEGER IWORK( * )
-
- REAL A( LDA, * ), B( LDB, * ), S( * ), WORK( * )
-
- IIIIMMMMPPPPLLLLEEEEMMMMEEEENNNNTTTTAAAATTTTIIIIOOOONNNN
- These routines are part of the SCSL Scientific Library and can be loaded
- using either the -lscs or the -lscs_mp option. The -lscs_mp option
- directs the linker to use the multi-processor version of the library.
-
- When linking to SCSL with -lscs or -lscs_mp, the default integer size is
- 4 bytes (32 bits). Another version of SCSL is available in which integers
- are 8 bytes (64 bits). This version allows the user access to larger
- memory sizes and helps when porting legacy Cray codes. It can be loaded
- by using the -lscs_i8 option or the -lscs_i8_mp option. A program may use
- only one of the two versions; 4-byte integer and 8-byte integer library
- calls cannot be mixed.
-
- PPPPUUUURRRRPPPPOOOOSSSSEEEE
- SGELSD computes the minimum-norm solution to a real linear least squares
- problem: minimize 2-norm(| b - A*x |)
- using the singular value decomposition (SVD) of A. A is an M-by-N matrix
- which may be rank-deficient.
-
- Several right hand side vectors b and solution vectors x can be handled
- in a single call; they are stored as the columns of the M-by-NRHS right
- hand side matrix B and the N-by-NRHS solution matrix X.
-
- The problem is solved in three steps:
- (1) Reduce the coefficient matrix A to bidiagonal form with
- Householder transformations, reducing the original problem
- into a "bidiagonal least squares problem" (BLS)
- (2) Solve the BLS using a divide and conquer approach.
- (3) Apply back all the Householder tranformations to solve
- the original least squares problem.
-
- The effective rank of A is determined by treating as zero those singular
- values which are less than RCOND times the largest singular value.
-
- The divide and conquer algorithm makes very mild assumptions about
- floating point arithmetic. It will work on machines with a guard digit in
-
-
-
- PPPPaaaaggggeeee 1111
-
-
-
-
-
-
- SSSSGGGGEEEELLLLSSSSDDDD((((3333SSSS)))) SSSSGGGGEEEELLLLSSSSDDDD((((3333SSSS))))
-
-
-
- add/subtract, or on those binary machines without guard digits which
- subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or Cray-2. It could
- conceivably fail on hexadecimal or decimal machines without guard digits,
- but we know of none.
-
-
- AAAARRRRGGGGUUUUMMMMEEEENNNNTTTTSSSS
- M (input) INTEGER
- The number of rows of A. M >= 0.
-
- N (input) INTEGER
- The number of columns of A. N >= 0.
-
- NRHS (input) INTEGER
- The number of right hand sides, i.e., the number of columns of
- the matrices B and X. NRHS >= 0.
-
- A (input) REAL array, dimension (LDA,N)
- On entry, the M-by-N matrix A. On exit, A has been destroyed.
-
- LDA (input) INTEGER
- The leading dimension of the array A. LDA >= max(1,M).
-
- B (input/output) REAL array, dimension (LDB,NRHS)
- On entry, the M-by-NRHS right hand side matrix B. On exit, B is
- overwritten by the N-by-NRHS solution matrix X. If m >= n and
- RANK = n, the residual sum-of-squares for the solution in the i-
- th column is given by the sum of squares of elements n+1:m in
- that column.
-
- LDB (input) INTEGER
- The leading dimension of the array B. LDB >= max(1,max(M,N)).
-
- S (output) REAL array, dimension (min(M,N))
- The singular values of A in decreasing order. The condition
- number of A in the 2-norm = S(1)/S(min(m,n)).
-
- RCOND (input) REAL
- RCOND is used to determine the effective rank of A. Singular
- values S(i) <= RCOND*S(1) are treated as zero. If RCOND < 0,
- machine precision is used instead.
-
- RANK (output) INTEGER
- The effective rank of A, i.e., the number of singular values
- which are greater than RCOND*S(1).
-
- WORK (workspace/output) REAL array, dimension (LWORK)
- On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
-
- LWORK (input) INTEGER
- The dimension of the array WORK. LWORK must be at least 1. The
- exact minimum amount of workspace needed depends on M, N and
-
-
-
- PPPPaaaaggggeeee 2222
-
-
-
-
-
-
- SSSSGGGGEEEELLLLSSSSDDDD((((3333SSSS)))) SSSSGGGGEEEELLLLSSSSDDDD((((3333SSSS))))
-
-
-
- NRHS. As long as LWORK is at least 12*N + 2*N*SMLSIZ + 8*N*NLVL +
- N*NRHS + (SMLSIZ+1)**2, if M is greater than or equal to N or
- 12*M + 2*M*SMLSIZ + 8*M*NLVL + M*NRHS + (SMLSIZ+1)**2, if M is
- less than N, the code will execute correctly. SMLSIZ is returned
- by ILAENV and is equal to the maximum size of the subproblems at
- the bottom of the computation tree (usually about 25), and NLVL =
- MAX( 0, INT( LOG_2( MIN( M,N )/(SMLSIZ+1) ) ) + 1 ) For good
- performance, LWORK should generally be larger.
-
- If LWORK = -1, then a workspace query is assumed; the routine
- only calculates the optimal size of the WORK array, returns this
- value as the first entry of the WORK array, and no error message
- related to LWORK is issued by XERBLA.
-
- IWORK (workspace) INTEGER array, dimension (LIWORK)
- LIWORK >= 3 * MINMN * NLVL + 11 * MINMN, where MINMN = MIN( M,N
- ).
-
- INFO (output) INTEGER
- = 0: successful exit
- < 0: if INFO = -i, the i-th argument had an illegal value.
- > 0: the algorithm for computing the SVD failed to converge; if
- INFO = i, i off-diagonal elements of an intermediate bidiagonal
- form did not converge to zero.
-
- FFFFUUUURRRRTTTTHHHHEEEERRRR DDDDEEEETTTTAAAAIIIILLLLSSSS
- Based on contributions by
- Ming Gu and Ren-Cang Li, Computer Science Division, University of
- California at Berkeley, USA
- Osni Marques, LBNL/NERSC, USA
-
-
- SSSSEEEEEEEE AAAALLLLSSSSOOOO
- INTRO_LAPACK(3S), INTRO_SCSL(3S)
-
- This man page is available only online.
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
- PPPPaaaaggggeeee 3333
-
-
-
-